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Abstract: The matrix model of ABJM theory can be formulated in terms of a Fermi gas in an 
t-h external potential. We show that, in this formalism, vevs of Wilson loops correspond to averages 

of operators in the statistical-mechanical problem. This makes it possible to calculate these 
vevs at all orders in 1/N, up to exponentially small corrections, and for arbitrary Chern-Simons 
coupling, by using the WKB expansion. We present explicit results for the vevs of 1/6 and the 
1/2 BPS Wilson loops, at any winding number, in terms of Airy functions. Our expressions are 



shown to reproduce the low genus results obtained previously in the 't Hooft expansion. 
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1. Introduction 

Localization techniques in superconformal field theories have provided matrix model representa- 
tions for partition functions and Wilson loop vacuum expectation values (vevs) on spheres. For 
M = 4 super Yang-Mills theories, these techniques were developed in [51], providing a proof of 
previous conjectures in [19, 14] which proposed a Gaussian matrix model formula for the vev 
of a 1/2 BPS, circular Wilson loop. This was extended to Chern-Simons-matter theories in 
[37, 35, 28]. In particular, a matrix model was obtained in [37] which calculates the partition 
function and the vev of the 1/6 BPS Wilson loops for ABJM theory [2] constructed in [17, 11, 52]. 
1/2 BPS Wilson loops were constructed and localized in [17], and their vevs are calculated by 
computing the averages of supertraces in the ABJM matrix model of [37]. 

Once the matrix models have written down, an important question is to extract from them 
the 1/N expansion of the observables, in order to test predictions based on the AdS/CFT corre- 
spondence. In the case of the 1/2 BPS Wilson loop of N = 4 super Yang-Mills, this is relatively 
straightforward, since the matrix model is a Gaussian one. In particular, in [14], a procedure was 
presented to obtain the full 1/N expansion of the 1/2 BPS Wilson loop, and explicit expressions 
were obtained for the leading term in the 't Hooft parameter, at all orders in 1/N. This term 
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gives, in the AdS dual, the leading contribution coming from strings with one boundary and 
arbitrary genus. 

The ABJM matrix model is much more complicated than the Gaussian matrix model. How- 
ever, its free energy can be computed to any desired order in the 't Hooft 1/N expansion [15, 16], 
in a recursive way. This is achieved by using the holomorphic anomaly equations of topological 
string theory [7], as adapted to matrix models and local geometries in for example [33, 20, 26]. 
For Wilson loops, results at low genus can also be obtained from matrix model techniques. The 
exact planar result was obtained in [45], and the first 1/N correction was calculated in [15] by 
using the results of [4]. In principle, one can compute higher genus corrections by using for 
example the topological recursion of [21], but this procedure becomes rapidly quite cumbersome. 
Unfortunately, we lack an efficient holomorphic anomaly equation for open string amplitudes 
which makes possible to go beyond the very first genera. 

In the context of ABJM theory, understanding the full 1/N expansion is however of great 
interest, since this gives quantitative information about the M-theory AdS dual. Equivalently, 
one can try to compute the observables in the so-called M-theory expansion. In this expansion, 
one still considers the limit of large N but k (the Chern-Simons coupling, or equivalently the 
inverse string coupling constant) is fixed. In [23], building on the results of [15, 16], it was shown 
that the full 1 /N expansion of the partition function could be summed up into an Airy function, 
after neglecting exponentially small corrections. This raises the question of finding a method 
for analyzing the matrix model directly in the M-theory regime, without having to resum the 't 
Hooft expansion. The method developed in [31] works directly in the M-theory regime and can 
be applied to a large class of Chern-Simons-matter theories, but in its current form it is only 
valid in the strict large N limit. 

A systematic method to analyze the matrix models arising in N > 3 Chern-Simons-matter 
theories, in the M-theory expansion, was introduced in [46]. The basic idea of the method is to 
reformulate the matrix model partition function, as the partition function of a non-interacting, 
one-dimensional Fermi gas of N particles, but with a non-trivial quantum Hamiltonian. In this 
reformulation, the Chern-Simons coupling k becomes the Planck constant h, and the M-theory 
expansion corresponds to the thermodynamic limit of the quantum gas. It was shown in [46] that 
the partition function of the gas could be computed, at all orders in 1/N, by doing the WKB 
approximation to next-to-leading order (neglecting exponentially small corrections). This makes 
it possible to re-derive the Airy function behavior found in [23] for ABJM theory, and generalize 
it to a large class of N = 3 Chern-Simons-matter theories. The Fermi gas approach provides as 
well an elementary and physically appealing explanation of the famous TV 3 / 2 scaling predicted in 
[39] and first proved in [15]: it is the expected scaling for a Fermi gas with a linear dispersion 
relation and a linear confining potential. 

In this paper we extend the Fermi gas approach of [46] to the calculation of vevs of 1 /6 and 1/2 
BPS Wilson loops. As expected, the vevs correspond, in the statistical-mechanical formulation, 
to averages of n-body operators. Since the gas is non-interacting, this can be reduced to a 
quantum-mechanical computation in the one-body problem, which can be in principle done in 
the semiclassical expansion. However, in this case a precise determination of the vev requires the 
resummation of an infinite number of quantum corrections. This is not completely unexpected: 
already in the calculation of the partition function in [46], there is an overall factor which is a 
non-trivial function of k and involves a difficult, all-order calculation of quantum corrections. For 
Wilson loops in ABJM theory, one can actually perform the resummation directly, and obtain 
a closed formula for the 1/N expansion of the 1/6 and 1/2 BPS Wilson loops in terms of Airy 
functions. In the case of a 1/2 BPS Wilson loop in the fundamental representation, the result 
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for the normalized vev is particularly nice, 



r i/2. I /2n^ Ai 



c-v*(n- 



A _ A 

24 3k 



W > = 4 CSC ( T ' r_ .J : 7TT (1-1) 



Ai 



f'- l/S (*-fi-& 



where 



C =Z2T- (1-2) 



2 

7T 2 /c" 

This result is exact at all orders in 1/iV, up to exponentially small corrections in N (corre- 
sponding to world sheet or membrane instanton corrections). The denominator in (1.1) is the 
partition function of ABJM theory as computed in [23, 46]. The corresponding expression for 
the 1/6 BPS Wilson loop, and for arbitrary winding, is more involved, and it is given below in 
section 4 (eq. (4.110)). 

The paper is organized as follows. In section 2, we start with a brief review of certain aspects 
of ABJM matrix model, and in particular, we review and extend the results of matrix model 
computations for the 1/6 and 1/2 BPS Wilson loop expectation values at genus zero and one. In 
section 3, we first proceed by recalling some standard techniques of quantum statistical mechanics 
in phase space, which are going to be used later on in this paper. We then turn into a brief 
review of the Fermi gas approach which was originally introduced in [46]. Section 4 is the core 
of our paper. We first demonstrate in 4.1 how we can include the Wilson loops in the Fermi gas 
formalism. We continue in section 4.2 by first computing the full quantum corrected Hamiltonian 
of the fermionic system, and then by calculating the corresponding Wigner-Kirkwood corrections 
for the quantum mechanical averages. In 4.3, we deal with the integration over the quantum 
corrected Fermi surface, and section 4.4 contains the explicit results for Wilson loop vevs and 
a detailed comparison with the 't Hooft expansion in the strong coupling regime. Section 5 is 
devoted to conclusions and prospects for future work. In Appendix A, we present the details 
of the matrix model computation for the 1/6 BPS Wilson loop correlator at arbitrary winding. 
Appendix B summarizes the results of the 'tHooft expansion at genus three and genus four. 

2. Wilson loops in ABJM theory 
2.1 1/6 BPS and 1/2 BPS Wilson loops 

The ABJM theory [2, 5] is a quiver Chcrn-Simons-matter theory in three dimensions with gauge 
group U(N)j. x U(N)_k and N = 6 supersymmetry. The Chern-Simons actions have couplings 
k and —k, respectively, and the theory contains four bosonic fields Cj, I = 1, ■ ■ ■ ,4, in the 
bifundamental representation of the gauge group. One can construct an extension of this theory 
[3] with a more general gauge group U(N\)k x U(N2)-k, but we will not consider it in detail in 
this paper. The 't Hooft parameter of this theory is 

N 

A=-. (2.1) 

A family of Wilson loops in this theory has been constructed in [17, 11, 52], with the structure 

Wj/ 6 = Tr R V exp J (iA^ + y WMfacA ds , (2.2) 

where is the gauge connection in the U(N)k gauge group of the first node, x fl (s) is the 
parametrization of the loop, and Mj is a matrix determined by supersymmetry. It can be chosen 
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so that, if the geometry of the loop is a line or a circle, four real supercharges are preserved. 
Therefore, we will call (2.2) the 1/6 BPS Wilson loop. A similar construction exists for a loop 
based on the other gauge group, and one obtains a Wilson loop associated to the second node 

W^ 6 . (2.3) 

In [37] it was shown, through a beautiful application of localization techniques, that both 
the vev of (2.2) and the partition function on the three-sphere can be computed by a matrix 
model (see [44] for a pedagogical review) . This matrix model is defined by the partition function 

Z^u(N,g s ) = — 2 Ultw „ ,*L-uA ~ e 29s(ElMl E " } ' ^ 
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where the coupling g s is related to the Chern-Simons coupling k of ABJM theory as 

9s = f- (2.5) 

One of the main results of [37] is that the normalized vev of the 1/6 BPS Wilson loop (2.2) is 
given by a normalized correlator in the matrix model (2.4), 

K /6 ) = (TV R (e^)) ABJM , (2.6) 

Notice that the Wilson loop for the other gauge group, 

(W 1 R /6 ) = (Tr R (e«)) AB3M (2.7) 

can be obtained from (2.6) simply by conjugation, or equivalently, by changing the sign of the 
coupling constant g s — > —g s - From now on we will then focus, without loss of generality, on the 
Wilson loop associated to the first node, and we will also assume that k > in the first node. 

The Wilson loop (2.2) breaks the symmetry between the two gauge groups. A class of 1/2 
BPS Wilson loops was constructed in [18] which treats the two gauge groups in a more symmetric 
way (see also [43]). These loops have a natural supergroup structure in which the quiver gauge 
group U(N) x U (N) is promoted to U (N\N), and they can be defined in any super- representation 

1/2 

1Z. In [18] it has been argued that this 1/2 BPS loop, which we will denote by , localizes to 
the matrix model correlator 

(<^< St -(t-e.)L M (2 - 8 » 

in the ABJM matrix model. Here, Str^ denotes a super-trace in the super-representation 1Z. In 
order to write this in more down-to-earth terms, we note that a representation of U(2N) induces 
a super-representation of U(N\N), defined by the same Young tableau 1Z (see for example [6]). 
Therefore, (2.8) can be also written as [6] 

-*(T_^?^f n 

k £ ' (2.9) 
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In this equation, which is the supergroup generalization of Frobenius formula, k = (kg) is a vector 
of non-negative, integer entries, which can be regarded as a conjugacy class of the symmetric 
group, xn{k) is the character of this conjugacy class in the representation 1Z, and 
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We will be particularly interested on Wilson loops with winding number n, which in the basis of 
representations are defined by 

n-l 
s=0 

Here, R njS is a "hook" representation with n boxes in total, n — s boxes in the first row, and 
one box in the remaining rows. For n = 1 we recover the usual Wilson loop in the fundamental 
representation. In terms of matrix model vevs, 

(W^) = (Tr(e n ^)) AB3M . (2-12) 

In view of (2.9), the 1/2 BPS Wilson loop with winding n is simply given by 

(W}J 2 ) = (W B V«> - (-l) n (W^). (2.13) 

In general, as it is clear from (2.9), the vevs of 1/2 BPS Wilson loops can be obtained if we know 
the vevs of 1/6 BPS Wilson loops, but the former are much simpler. 

2.2 The geometry of the ABJM theory 

In [45] the correlation functions of the ABJM theory that are accessible via localization [37], 
i.e. the partition functions and the Wilson loops, are mapped via the spectral curve of the 
Chern-Simons matrix model, to geometric invariants of the elliptic curve 



1 

<fiX • ' • <f%Y 



H(X,Y) = X + ^^+Y + -^ + l = (2.14) 



with are related to meromorphic differentials of the third kind. In particular in the planar 
limit the partition function and the Wilson loop expectation values of ABJM are simply related 
to closed periods of these differentials. The higher N corrections are related by a recursive 
procedure, which amounts to integration of the loop equation, to these periods. In (2.14) X, Y 
are C* variables and eqn. (2.14) is the B-model mirror curve of the local Calabi-Yau geometry 
M cy = 0(—K F i xF i) — > P 1 x P 1 , i.e. the total space of the anti canonical line bundel over P 1 x P 1 . 

After multiplying (2.14) with XY, homogenizing it to a cubic with W, swaping W with — Y 
and rescaling X i-> X (p\ one gets the curve 

H(X, Y) = Y 2 - Y (l + A>i + X 2 ^j + X 2 = (2.15) 

One might parameterize the C* variables X = e u and Y = e v . Then the relevant meromorphic 
differentials of the third kind are given 

f i k = ve ku du = \og(Y)X ( - k -^dX, A; = 0,1,..., (2.16) 
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where 

Y= a{ ^± 1 2 ^iX)- (2-17) 

This form is typically for local mirror geometries and with the above parameterization the dis- 
criminant is given as 

4 2 

<r(X) = Y[(X - Xi ) = a{X) 2 - AX 2 , with a(X) = (1 + Xtpx + X 2 ^) . (2.18) 
i=l ^2 

The branch points involve square roots of the tpi, but with an appropriate ordering one has 

4 

11 1 1 

(pi = -~ y^Xi, — = -(xi + x 2 -x 3 -x 4 ), xi = — =:a, x 4 = — =: -6 . (2.19) 

2 <P2 4 x 2 x 3 

i=i 

Note that (2.15,2.17) defines the same family of (hyper) elliptic curves as 

y 2 = a(x) , (2.20) 

where we identified X, Y with x, y. This identification amounts to a compactification of C* 
variables X, Y and does not affect integrals over closed cycles up to one important subtlety. At 
X — > oo, no behaves like 

In the compactifications one has to regularized the form to Ho( x ) = [io(X)\ x= x — |log(x). 
Derivatives of the /j.q(x) w.r.t. to tpi are related to standard elliptic integrals on (2.20). 

The Xi = ^ parameters of the Chern-Simons matrix model are mirror coordinates and as 
such identified with the periods 

1 



A; = — / HO- (2-22) 
J c . 

over the C-cycles the geometry 

C 1 = (l/a,a), C 2 = (-6,-1/6) . (2.23) 
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The homology relations imply that the Cj periods are non identical because of the pole in fi n . 
In particular for fj,Q it is clear from Figure 1 and (2.21) that there is an exact relation between 
the periods (2.22) 

exp(27ri(Ai - A 2 )) = — • (2.24) 

For this reason the ABJM slice Ai = A2 mod Z can be identified with an algebraic sub manifold 
of complex deformation space of (2.14) given simply by 

fi = ¥2 = <P = in. (2.25) 

In particular on the slice one has 

d^Hk = w fe = dx , (2.26) 

y/<r{x) 

i.e. all closed integrals of fik on (2.15,2.17) are determined up to a constant by standard elliptic 
integrals on (2.20). For latter reference we note that the parameterization of the branch points 
by k is 

a(/e) = - (2 + m + V^(4i - «)) , 6(k) = ^ (2 - Ik + v/-K(4i + «)) . (2.27) 

On the slice (2.15) is an algebraic family of elliptic curves r°(4) curves with j-invariant 

_ 16 - 16p 2 + ^ 
3 ~ 1728^(16 -v? 2 ) ' 1 8j 

The curve is related to the Tq(2) curve of pure SU (2) SW-theory y 2 = (x 2 — u) 2 —A 4 by identifying 
u = ±(1 — ^-)A 2 . Indeed the periods integrals of (j,q are annihilated by a single Picard-Fuchs 
differential operator of topological string theory for M cy identifying both Kahler classes of the 
P 1 . i.e. T 1 = T 2 . It reads 1 

V = {^ 2 e% - 16(2^ - l) 2 )6^ = vV(^ + l) 2 + 160j)0 v = tpVwdy, (2.29) 

where X = x-^ is the logarithmic derivative. T>hoi annihilates the periods over the holomorphic 
differential wo on (2.20) as a consequence of (2.26). This differential equation has (2.29) three 
critical points (p 2 = 0, ip 2 = 16 and (p 2 = 00. Let us describe the behavior of the periods at these 
points and determine the analytic continuations and the monodromy action. The weak coupling 
point of ABJM is the point (p = 0. In the w = p 2 variable the period basis looks like 




n = L an = d\F , with 





/ 1 \ 

. (2.30) 



& + ^- 2 +0{w 2 )] 



where d\F° = F\ — (2 + ib)\ — \ and b = 41og ^ 2 ) +1 , The recursion of A can be summed to yield 

A= 8^ 3jF2 G'^ ;1 4 ; "T6) • (2 - 31) 



lr rhe formulas 8 X = a6 y if y = x a and [8 X , x z ] — ax a make the comparison with [41] trivial. 
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This function plays the role of the mirror map at the orbifold while d\F® is the dual period. The 
pair defines by special geometry F® and the polarization on the ABJM slice. 

The point (p 2 = oo is the strong coupling point of ABJM A — > oo and corresponds to the 
large radius point of topological string. The topological string basis is obtained by the local limit 
of a compact Calabi-Yau manifold and is half integral in the homology of the curve (2.15) 

/ 1 0\ 

n= -l i o n ts . (2.32) 
V o -I i/ 

In the z = ^2 coordinates the topological string basis reads 

n ts =f T ) = ( ^ l [log(z)+Az + 0(z 2 )} ]. (2.33) 

\9tF° w J \-$ (2^) * Pog'CO + 4^1og(z) + 8z + 0{z 2 )] -±J 

Here Q = exp(2mT) and d T F° w can be integrated to the genus prepotential Fg W (Q) = —qT 3 — 
j^T + c + YldLi n dLi3(<? d ), the generating function of g = BPS invariants summing over the 
degrees d\ + d2 = d w.r.t. to both Kahler classes of the P^s 2 . The same basis at the conifold in 
the u = (1 — ^) coordinates 3 



1 



n 



ts 



d TcF° c = ^ Il [^T c \og{u) + lu 2 + 0{u 3 )] + 2biT c + c . (2.34) 



From this we get the T°(4) monodromies in the II basis 

/ 1 0\ /l 0\ / 1 

M v=0 = -1 -1 -4 , M^2 = _ 16 = 3 4 , M^ =oc = 1 1 | . (2.35) 

\ o o -l / \o-i-i/ \-i 

One checks (M ¥ ,2 = _ 16 M (/3=00 )~ 1 = M^ =0 . 

In topological string theory or N=2 4d supersymmetric gauge theory, the coupling constants 
are complex. At the various critical points one has to chose appropriate coordinates, which are 
either invariant or reflect invariances of the theory under the local monodromy. For example at 
large radius or the asymptotic free region of the gauge theory one can chose T as variables and 
the monodromy T — > T + 1 is understood as a shift in the NS B field of topological string or the 
(9-angle of Yang-Mills theory. The canonical choices of other coordinates in different regions in 
the moduli space correspond to a change of polarization. 

Because in ABJM theory the coupling constant is real there is a priori no need to consider 
the action of the monodromy. The polarization is picked once and for all at the weak coupling 
point. The choice made here is identical to the one made in topological string theory at this 
point in moduli space. 

However this polarization is not the one of the topological string at large radius. The coupling 
of the ABJM theory A behaves at large radius like 

A = d T F®, - l T = --T 2 - -T - — + 0(Q) . (2.36) 
2 2 2 12 v ; K J 



Up to a constant c = 2 ^^ 3 , which depends on the regularized Euler number of the local geometry for x — 4. 
3 The irrational constant c = .3712268727. . . is fastest iterated using the Meijers function [41]. 
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It is crucical to obtain the famous iVz scaling of the free energy energy to integrate the 
special geometry relation from £?-cycle integral d\F^ with respect to A [15]. This yields 4 



9s 



fc 2 Ai + O(X ,e 



■2w\ 



(2.37) 



The relation of the topological string theory to the AB JM theory at this point is therefore given 
by a change of polarization. 

What is remarkable is that despite the fact that the action of the monodromy is not under- 
stood in ABJM theory, the higher genus contributions to the partition function of the theory 
have the same modular invariance under T°(4) as the one of the topological string. One might 
speculate that the monodromy at the strong coupling region reflects an, so far not understood, 
invariance of the theory, which involves non-perturbative effects. Note that this monodromy does 
not change the leading Nz behavior. A related topic is the g-Wilson loop itself. It corresponds 
to an integral over the C cycle, while integral of the same differential over the dual B-cycle has 
no interpretation in ABJM theory. Again if the monodromy has a meaning ABJM it would mix 
the two types of cycles. 

2.3 Wilson loops in the geometric description 

The exact planar vevs of 1/2 BPS and 1/6 BPS Wilson loops were obtained in [45], from the 
exact solution of the ABJM matrix model at large N. We will now review (and extend) these 
results. 

The planar limit of the matrix model is completely determined by the densities of eigenvalues 
in the cuts, which were also obtained explicitly in [45]: 



Pi(X)dX 
P2(Y)dY 



1 



2vr 2 Az 
1 

2n 2 M 



tan 



tan 



laX - 1 - X 2 
/3X + 1 + X 2 



aY - 1 - Y 2 



dX 

~X' 

dY 
~Y' 



(2.38) 



where 



a = a H — , (3 = b + -. 
a b 



(2.39) 



These densities are normalized in such a way that their integrals over the cuts are equal to one. 
We drop this normalization factor in the following formulas. It can be reinstalled after taking 
the derivatives. 

Keeping track about the residuum at X = oo analogously to (2.21) we can write simpler 
expressions for the densities, valid in the compactified variables x, y. The planar | Wilson loops 
read in terms of those 



(Wh 9 =0 = ^f »n, (Wh 9 =0 = ("IT ^2 



(2.40) 



id ™ Jc 2 

while the planar \ Wilson loops follows by the relation in homology as the 7-period, i.e. the 
residuum at infinity 

(wh g =o = ^fpn. (2.41) 



4 As further explained in [15] it is natural to shift A by A — y A 
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Since the u: n = du>[/, n , n = 0, . . . , are not independent elements of the cohomology of the curve, one 
can related all Wilson-Loop expectation values to the integrals of Let us denote R n (^p) = J y 
the residue of fj, n at x = oo then we get a relation in homology involving the 

£ n u) -oo n = d ip R n (ip)x 3 dx , where C n = p^^d^ +p°(<p) • (2.42) 

The coefficients p° (<p) and p^if) are polynomial in 92 and given by the Griffith reduction method. 
For the first few we get 



p°l 


— <P 
~ 4 ' 






Pi 


v\ 


= 1, 






P\ 


P3 


9<fi 
4 






P\ 


Pa 


= 1 + 


i(V 

3 


5<p 4 
24 ' 


P\ 



0, Ri = \ip, 



(2.43) 



1 

This relates (Wn) g =o to A, e.g. 



{Wh 9 =* = \ J K\{n)dn+ l -n. (2.44) 

The integration constant is zero as fj,i has no constant residuum. 

The relations (2.42) are homological relations, which imply a differential relation of the B- 
cycles integrals over jx n to d\F°. Since A, d\F° are related by special geometry (2.42) imply for 
each k differential relations between the two kind of Wilson-Loop integrals, which can be view 
as an extension of special geometry to the Wilson loop integrals. 

The Wilson loop vevs have a genus expansion of the form, 

00 

« /6 ' 1/2 ) = E^ 1 « /6 ' 1/2 )^ (2-45) 

and this expansion maps to the genus expansion of the ABJM matrix model correlators (2.12). It 
is a standard result in matrix model theory that the planar correlators of this type are moments 
of the eigenvalue densities, 

(Tre^% =0 = [ P i(X)X n dX. (2.46) 
JCi 

In this form, these correlators are difficult to compute, but as in [9], their derivatives w.r.t. k 
are easier to calculate and given by 

where 

1„ = f X " dX (2.48) 

2 J J(aX-l-X 2 )(J3X + l + X 2 ) 

1/a 

which can be calculated in terms of elliptic integrals. The computation for n = 1 was done in 
[45], and in Appendix A we compute them for a positive integer n. In order to make contact 
with the Fermi gas approach, where subleading exponential corrections are neglected, we want 
to extract their leading exponential behavior. One finds, 

"ft fl — X / ' \ 

2n~ HogK-y-fln-iJ (2.49) 
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where 

Hn = jz\ ( 2 - 5 °) 



d 

d=l 



are harmonic numbers (for n = 1, we set i?o = 0). It then follows that 

(*«-?- "-OHO?))- < 2 - 5i » 

From this we deduce that, for the 1/2 BPS Wilson loop, one has 

(H^_-!^( 1+0 ^)). ,, 52 ) 

This agrees with a result obtained in section 8.2 of [15], where the generating function of these 
vevs, with an extra 1/n factor, was shown to be a dilogarithm in the variable in. 

The regime of large n corresponds to the regime of large 't Hooft coupling [45] , and one has 
from (2.31), 

W-^ + ^+of-Vi, (2-53) 



2vr 2 24 \k 2 , 
which is immediately inverted to 

k = e"^ (l + O (V 27rv/ ^) ) . (2.54) 

It follows that the 1/6 Wilson loops go like, 

(W„ 1/6 ) » e n7rV ^, (2.55) 

and for n = 1 this is in agreement with the AdS calculation in terms of fundamental strings 
[17, 11, 52]. 

2.4 The higher genus calculation from the spectral curve 

The essential information of the higher genus expansion is encoded in the expansion of the 
resolvent 

oo 

M{x) = Y J 9l 9 $\x) • (2-56) 

9=0 

The densities for the Wilson line integrals of winding k at genus g are then 

$\x) = x k $\x). (2.57) 

The main task is hence to determine the expansion (2.56) by the matrix model recursion. The 
simplest formulation of the latter uses the hyperelliptic curves 

V 2 = <t(x) (2.58) 

and as meromorphic differential defining the filling fractions 

$(p) = y(p)dx(p) . (2.59) 
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Instead of (2.59) we want to work with /io as differential denning the filling fractions. Let us 
denote the points on the ±-branch of (2.20) p and p. I.e. both points map to the same x value. 
In the recursive formalism the discontinuity of <3?(p) — <fr(p) over the cuts is essential. Likewise at 
a(X) = the two branches of the curve (2.17) come together. The difference of fio on the two 
branches is however 

4*(P) - = ^rtaoh- (igp) ix . (2 . 60 ) 

One can now redefine y in order to match these differences. This leads to the definition of a 
curve 

f = M 2 (x)a(x) , (2.61) 
on which (2.59) is equivalent to [iq on (2.15,2.17). Note that moment function 



M{x) = l tanh- 1 ^^ (2.62) 

x^{x) a(x) 

does not modify the branch points. In particular it does not introduce new ones. The moment 
function appears in the recursive formulation of Akemann [4] and can also be implemented in 
the Eynard Oratin formalism [21]. 
The recursion formula reads 5 

d c l( ) f g \ 
Wg(p,p!,...,p k ) = ^2 Res q=Xi - ' ^2 ^2 W h (q,pj)W g - h (q,p K \j) +W g - 1 (q,q,p K ) . 

i y { q > \h=0JcK J 

(2.63) 

Here K, J are index sets K = {1, . . . , k} etc. In principle we are only interested in the W g (p ), 
however the higher genus recursion requires for g = 2 to calculate amplitudes with up to three 
legs at genus 0. The dS(p, q) are the unique meromorphic differentials with only simple poles at 
q = p and q = p, whose integral w.r.t. to q over the A cycles, which we call in our context C, 
vanish 

dS(q,p) 9 ~ p dS(q,p) q ~ p - * f dS(q,p) = . (2.64) 

y P IP JqCAi 

The decisive technical tools to solve the recursion are the so called kernel differentials 

x <^ Res „(^_J_). 

Multiplying an expression f(q,Pi,Xi) by dS y^ and taking the sum of the residua at q = Xi is 
the crucial step in solving the recursion, so let us denote 

@(p,q)f{q,Pi,Xi) = Y^ Res i=^ { ^ S ~^ f(viPii x i^ ■ ( 2 - 66 ) 

Beside the genus zero resolvent Ho(x) to start the recursion one needs the annulus amplitude 
W (p,q) = - , 1 + ^_ ^EL ^ + *P> q }_ . (2.67) 



5 One writes W g (pi, . . . ,p k )dpi ...dp k := ui g (pi, . . . ,p k ) = H (pi, ■ ■ ■ ,Pk 
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which is related to the Bergman kernel by 

B(p, q) = (w (p, q) + ^^^2 ) d P d( l ■ ( 2 - 68 ) 

On an elliptic surface A(p, q) as well as the kernel differentials are given in terms of elliptic 
integrals. 

A(p,q) = {p-xi){p- x 2 ) + {p- x 3 )(p- x A ) + (xi - x 2 ){x A - x 2 )G(k) (2.69) 

where 

h 2 = 

(Xl - X 3 )(x 2 ~ Xa) 



k2 = ( Xl - x 2 )(x 3 - x 4 ) 



is the elliptic modulus and 



G{k) = |^ (2.71) 



K{k) 

is the ratio between the two complete elliptic integrals 

— — n 

E{k) = [ 2 \A - k 2 sin 2 6d9, K(k) = [* - . (2.72) 

Jo Jo VI -k 2 1 



sin 



2 , 



As explained in [40] the ordering of the branch points here follows the one appropriate for 
local P 1 xP 1 , which is obtained from the one in [4] by the exchange 

x 2 o x A . (2.73) 

The expression for the kernel differentials follows from a Taylor expansion of 

dS(p, q) 1 



y(q) M(q)y r ^p) \p 

around the branch points. Here 



7T^+^ (1) (g)) d P ( 2 - 74 ) 



= ,CC«( fl ) = -V^-^-^ cW{q) (2 . 75) 
is a normalization of the C (or equivalently the A) cycle integral 

C«(g) = / ~ d \ y= , (2.76) 
W Jc 2m (q - x)Jd V ' 

so that so that the last property (2.64) hold. Note that this integral has to be regularization if 
q approaches the branch points of the cuts defining the C cycle integral 



' J_ f dx 

2m JC (q-x)y^ 



if Xi is not a branch point defining C 

q=Xl I ( 2 - 77 ) 
2ln Ic J^Wa ~ if x * is a branch P° int defining C 



This definition of the regularization ensures that one can move the contour from the x\ — x 2 
cut to the X3 — X4 cut without getting a contribution from the poles. As a consequence the so 
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defined integrals C(xi) obey a symmetry under certain permutations of the branch points. We 
can evaluate e.g. the manifestly regular integral 6 



«4 =M W (X 4 ) 



(X 3 - Xj) 

(xi - x 4 ) 



G(k) + 1 



x 4 - x 3 

and obtain from the symmetrization the evaluation at the other branch points 

x 2<~* x 3 



(2.78) 



x 3 



^x 2 



Higher kernel differentials are therefore given by 

1 1 d n ~ l f 1 



X 



(n) 



— +M®{q) 



(2.79) 



(2.80) 



q=Xi 



(n-1)! yffiqjdtr- 1 [M(q) \p-q 

Here = /CC«( g ) and since the normalization factor K, is independent of q. the only 

nontrivial task is to calculate the dei 
do this. One fast way is to compute 

{-l) n ~ l {n- 1)! 



nontrivial task is to calculate the derivatives C^ n \q) = ^-i q C^(q). There are various ways to 



C^ n \q) 



2m 



dx 

C (q - x) n y/a 



(2.81) 



These integrals with poles at finite points and are very similar as the ones with poles at infinity 
and by similar formulas they can be expressed by linear expressions in K(k) and E(K) with 
rational coefficients in the moduli. In particular the normalized integrals J\f( n \q) depend only 
on the ratio of elliptic functions G(k) To get expressions which are valid at all branch points one 
calculates first M^ n \x 4 ), which is regular and uses (2.79) to get Af( n \xi). These derivatives have 
symmetric expressions in terms of the branch points and the a*. E.g. the first two derivatives 
are 



— Oti 



0C n X j 



15 



(Xj XiY 



{Xj Xj^^Xfc Xi 



Eventually one needs integrals over meromorphic forms with mixed poles 



(2.82) 



-.dx , 



(x — p) k \f a{x 
which are obtained from the obvious relations 

The genus one differential is then determined by evaluating 

W 1 (p) = e(p,q)W (q,q) 



(2.83) 



(2.84) 



(2.85) 



3 Here we make contact with the shorthand notation a, introduced in [4]. 
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using (2.67,2.65) as well as the explicite formulas for the kernel differentials for elliptic curves. 
With this one can order W\{P) according to it s poles in the branch points 



W±{p) 



E 



+ ^- + C i 



(2.86) 



where 



1 1 

16M( Xi )' 



Bi 



1 M'{xi) 1 



16M 2 (xi) 8M(xi) 
i 1 M'( Xi ) 



^— ' T - .- — 1"- / 



_ 1 1 v - Qj ~ Qj _ 1 M'^jJ , OLj 

1 ~ 48 M( Xi ) j£ Xj - Xi 16 M 2 ( Xi ) °* 8M fa) 



. X% Xj 
3+i 



(2.87) 



3/ ^ 37 7 



To obtain the form W2(f>) 

W 2 (p) = e(p,q) (W 1 (q,q) + W 1 (q)W l (q)) 

one needs W±(p,pi) from 

Wi(p,pi) = e(p,(z) (Wo(9,g,pi) + 2Wi( g )W ((z)) 

and Wo(p,p±,p 2 ) from 

Wi(p,pi,p 2 ) = 2@{p,q)W {q,Pi)W {q,p 2 ) . 



(2.88) 



(2.89) 



(2.90) 



By repeated application of the recursion one expresses any amplitude through a calculating of 
repeated residue of products of the annulus amplitude. E.g. 



W 2 (p) =20(p, q)@(q, qi)Q(qi, q 2 )W (q 2 , q)W (q 2 , 

20(p, 5)0( g , gi)0(9i, ©)Wb(gi, g)W («2, g 2 )+ 
6(p, q)e(q, qi)Q(q, 92)W (gi, qi)W (q 2 , q 2 ) . 



(2.91) 



It is easy to derive that for amplitude with genus g and h holes all terms will be of the general 
form 

Wg )h ~ e 2 »- 2+ft W^ h_1 . (2.92) 

However the number of terms grow exponentially with g and h. A few example for the number 
of contributions counted with multiplicity is given in the table below. 





9 





1 


2 


3 


4 


5 


h 
















1 




disk 


1 


5 


60 


1105 


27120 


2 




1 


4 


50 


960 


24310 




3 




2 


32 


700 


19200 






4 




12 


384 


12600 









Table 1: Number of terms involved in the recursive definition of W„ h- 
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Since W ~ G and each increases the power of G by one we get for the leading power 
Wg^h ~ Q3g+2h-3_ ]y[ ore precisely the W 9 ^ are meromorphic differentials with the following pole 
structure 

fZg-2+h h 4 ^0') \ 

W '^' • • • = wvm { g S § J • (2 ' 93) 

where the 

p=0 

are polynomials in the ratio of the complete elliptic integrals. The calculation at higher genus is 
obviously quite involved. Nevertheless one can make a general statement about the logarithmic 
structure of Wilson loop integrals at strong coupling. Since G(k) ~ \ g( K ) + we see 

that the highest inverse powers of jj^^y at leading order in k go with log (J)3 g -i • For the 1/6 
Wilson loop there will be a positive power of log(/c) at leading order due to the integration of 
the meromorphic differential 0J g {p) over the C cycle. The structure can be checked e.g. at genus 
two (4.122), 

For W2(p) we found an explicite expression for general moment functions. To write down all 
A*fl takes however several pages. We display the coefficient of the leading pole 

and af 9 2+h ^ := Y^i a ff multiplying the highest power of G, in the W2(p) case G 5 
(5) _ I5(x 1 - x 3 )(x 2 - x 4 ) 3 



2 



1 8M{x 1 )M{x 2 )M{x 3 )(x 1 - x 2 ) 4 (xi - x 4 ) 2 (x 2 - x 3 ) 4 (x 3 - x A ) m) 

3(xi -x 3 ) 3 (x 2 -x 4 ) 

+ 16(xi - x 2 f( Xl - x 4 ) 2 (x 2 - x 3 ) 4 M(xi)M 2 (x 2 ) + Sym • 

The other expressions are available on request. We obtain the weak coupling expansion of the 
1/2 Wilson line at genus 1 from (2.86) 

1. % 11. 97. 10331ivr 7 A 7 . . 

(W5)a=i = — wrA nr 3 \ 3 + ztt 5 A 5 + (A 9 ) . (2.97 

x 19 6 36 360 30240 v ; v 1 

Similarly we get from 

i 1407HvrA 12319ivr 3 A 3 693989iyr 5 A 5 59378093ivr 7 A 7 , 9 , 

^ 2 ' 9=2 ~ 120960 + 90720 14515200 + 152409600 + ^ '' ^' ' 

Note that (W2'6) s=lj2 are regular at the weak coupling point. This would be very hard to verify 
from the Fermi gas approach. 

On the other hand we can check that the strong coupling expansion of the {W^) g= i t2 Wilson 
loop matches the corresponding terms in the expression (4.119) (4.122). 
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3. The Fermi gas approach 



3.1 Quantum Statistical Mechanics in phase space 

The Fermi gas approach to the ABJM matrix models (and to other matrix models arising in 
N = 3 Chern-Simons-matter theories) is based on an exact equivalence with a quantum Fermi 
gas of N particles with Planck constant h = 2nk, and an evaluation of the different observables in 
the semiclassical expansion. For this reason, it is convenient to formulate the quantum mechanical 
problem in Wigner's formalism. We will now review some of the basic tools that we need to set 
up the formalism. 

We recall that, to construct the Hilbert space for a space of indistinguishable particles, one 
introduces a projection operator on totally (anti)symmetric states 

where 

7] = ±l (3.2) 
for bosons and fermions, respectively. This operator satisfies 

Pi = P v - (3.3) 

Let 

|Ai,---,Aiv) (3.4) 

be the basis of space eigenstates for an iV-particle system T~Ln of distinguishable particles. The 
appropriately (anti) symmetrized states 

|Ai, ■ - ■ ,\n} = VmPr,\\i,--- > X n) = -^=Y1 V <a) \K(i),--- ,K(N)) (3-5) 



constitute are a basis of the Hilbert space of bosons/fermions £>tv, Fn- The resolution of the 
identity in Bn/Fn reads 

^ J dA |Ai, ■ ■ • ,Xn}{\i,--- , Ajv| = 1. (3.6) 

A n-body operator O is an operator which is invariant under any permutation of the particles, 
and acts on a state oi%N as follows, 

0\\ 1 ---X N ) = ^ 0(A il ,---,A ifc )|Ai---A JV ) (3.7) 

• l<h^-^i k <N 
For example, for a one-body operator we simply have 

N 

0\\ 1 ---\ N )=J20(\ i )\\ 1 ---\ N ), (3.8) 

i=i 

where O(A) is an operator on the Hilbert space of a single particle. 
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In the canonical ensemble, the thermodynamic properties of the system are encoded in the 
canonical density matrix. For a system of distinguishable particles, the canonical density matrix 
is given by 

Pd({xi,--- ,x n }, {x[, ■ ■ ■ ,x' N };(3) = (xi ■ ■ ■ x N \e~^ H \x[ ■ ■ ■ x' N ), (3.9) 

where H is the total Hamiltonian of the N particles. For bosons (respectively, fermions), we 
have to (anti)symmetrize it in an appropriate way, to obtain [22] 

p({xi,--- ,x N },{x[,--- ,x' N };(3) = ^ V^PdHxi,--- ,x N },{x' a{1) ,- ■ ■ , x' a{N) }; p) 

= ^7 {xi ■ ■ ■ zjvle^lxi • • • x' N }. 

(3-10) 

In order to compute the vevs of many-body operators in the canonical ensemble, it is useful 
to introduce density submatrices or reduced density matrices (see for example [22, 57]). The 
reduced n-particle density matrix is defined as 

TV! f , , 

p n ({x!,--- ,x n },{x[,--- ,x' n };P) = ^ N J dx n+ i---dx N p({xi,--- , x N }, {x[, ■ ■ ■ ,x' N };/3) 

(3.11) 

The thermal average of an n-body operator O in the canonical ensemble is defined by 

(0) N = Tv(pO) (3.12) 

where we are using unnormalized vevs. This can be computed in terms of the n-reduced density 
matrix as follows 

(O)n = ~\ J dx 1 ---dx n O(x 1 ,-- ■ ,x n )p k ({x 1 , ■ ■ ■ ,x n }, {x\, ■ ■ ■ ,x' n };/3). (3.13) 

In our conventions, the canonical partition function is defined as the thermal average of the 
identiy, 

Z N = Tr(p). (3.14) 
We note that, in a system of non-interacting particles, the density matrix (3.9) factorizes, 

N 

Pd{{xi,--- ,x n },{x[,--- ,x' n };/3) = Y[p(x i ,x' i ), (3.15) 

i=i 

where p(x, x') is the canonical density matrix of the one-particle problem. 

In many situations it is more useful to work in the grand-canonical ensemble, where the 
reduced density matrix is defined as (see for example [57]) 

oo 

Pn°({xi,--- ,x n },{x'i,--- ,x' n };P,z) = z N Pn{{xi,--- , x n }, {x[, ■ ■ ■ ,x' n };/3) (3.16) 

N=n 

and as usual 

z = e^ (3.17) 
denotes the fugacity. The grand partition function is 

3 = 1+ ^ z N Z N , (3.18) 

JV=1 
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and the vev of an n-body operator in this ensemble can be simply expressed in terms of a sum 
of canonical vevs over all particle numbers, 



(O) 



GC 



N 



(3.19) 



N=n 



In the case of non- interacting gases, the grand-canonical density matrix has a very simple 
form (see for example [57, 30]): 



p£ C ({zi,--- ,*„},{*'!,■■■ ,x' n };p,z) = ~ Yl ^ (CT) n n (^ W ;M; ( 3 - 2 °) 



i=l 



where S is the grand-canonical partition function, and 



n(x, x'; P, z) = l^x 



z -i e /3H _ 



(3.21) 



is the occupation number operator in the position representation. The relationship (3.20) can 
be derived by using creation and annihilation operators [30]. There is also an elegant derivation 
in the case n = 1 by using the so-called Landsberg's recursion relation. This relation is based 
on the analysis of the sum over permutations in terms of conjugacy classes, and it was originally 
derived for the canonical partition function of ideal quantum gases (see for example [42]). It is 
however straightforward to generalize it to density matrices [12], and one finds 



v 



Pl (x,x';P) = J^r/ 1 p(x,x';£p)Z N _ e , 



(3.22) 



where p(x, x'; P) is the density matrix for the one-particle problem. We now sum over all -/V with 
the fugacity z N to obtain the grand-canonical, reduced density matrix, 



oo N 



pf c (x,x';p,z) = Y,Pi(^^';P)z N = V Y j 12 Zn -" z 



N-e 



N= 


i 






/ oo 


A 


E 






=0 




f 








z 



1 



x 



n=i e=i 

oo 



-epH 



x') fr? 



(3.23) 



z 1 e /3H — r] 

We conclude in particular that the vev of a one-body operator in the grand-canonical ensemble 
is given by 

O 



(O) 



GC _ - 



HTr 



(3.24) 



z -l e PH _ r?/ 

where the operator O appearing inside the trace is understood as the operator restricted to the 
one-particle Hilbert space. 

In order to calculate the quantum-mechanical averages, we will use a semiclassical or WKB 
expansion. The most convenient framework to do this is the phase space formulation of Quantum 
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Mechanics (see [32, 56] for detailed expositions). We first recall that the Wigner transform of an 
operator A is given by 



Aw(q,p) = dq' (q 



A 



q+ q -)j pq ' /h 



(3.25) 



The Wigner transform of a product is given by the ^-product of their Wigner transforms, 

(AB) w = A W *B W (3.26) 



where the star operator is given as usual by 

★ = exp 



ifr /tr 
~2 



and is invariant under linear canonical transformations. Another useful property is that 

dpdq 



TiA 



2irh 



A w (q,p). 



(3.27) 



(3.28) 



Let H be the Hamiltonian of a one-particle, one-dimensional quantum system, and let Hyj 
be its Wigner transform. Following [25] we notice that it is possible to expand any function f(H) 
of H around Hyj(q,p), which is a c- number. This gives, 



(3.29) 



r>0 



The semiclassical expansion of this object is obtained simply by evaluating its Wigner transform, 
and we obtain 

f(H) w = ^l i f (r) (Hw)Gr (3.30) 



r>0 



where 



w 



G r = [(H-H w (q,p)y 
and the Wigner transform is evaluated at the same point q,p. Of course, one has 

Go = 1, Gi = 0, 



(3.31) 



(3.32) 



and the quantities G r for r > 2 can be computed by using (3.26). They have an h expansion of 
the form 

Gr= Yj h2n ^ n) , r>2. (3.33) 

This means, in particular, that to any order in h 2 , only a finite number of Gr$ are involved. One 
finds, for the very first orders [25] , 



G2 

Gz 



4 
4 



d 2 H w d 2 H w _ f d 2 H w Y 
dq 2 dp 2 \ dqdp J 



( dH w \ 2 d 2 H w ( dH w \ 2 d 2 H w n dH w dH w d 2 H- 



\ dq ) dp 2 



+ 



dp 



w 



dq 2 



dq dp dqdp 



(3.34) 



0(h 4 
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One can then apply this method to compute the semiclassical expansion of any function of the 
Hamiltonian operator. A particularly important operator is the distribution operator 

n(E) = Q(E — H) , (3.35) 

where O(x) is the Heaviside step function. The trace of this operator gives the function n(E), 
counting the number of eigenstates whose energy is less than E: 

n{E) = Tvh(E) = J2@(E- E n ). (3.36) 

n 

One can regard also the operator (3.35) as the Fermi occupation number operator in the limit of 
zero temperature. When we apply (3.30) to (3.35), we find, 

00 (—Y] r 

h(E) w = Q(E - H w ) + ^-p-0 r a (r-1) (£ " H w ), (3.37) 

r=2 r ' 

and evaluating the trace according to (3.28) one obtains the useful formula, 

When (3.30) is applied to the canonical density matrix at inverse temperature (3, one finds, 

(•-"*) w - (f^ft) ( 3 - 39 ) 

We will call the the functions Q r appearing in (3.30) Wigner-Kirkwood corrections, and the 
resulting expansions Wigner-Kirkwood expansions. These corrections were originally introduced 
by Wigner and Kirkwood in their study of the semiclassical expansion (3.39) of the canonical 
partition function. Note that (3.39) can be interpreted as saying that the Wigner transform 
of the canonical density matrix is the generating function of the Wigner-Kirkwood corrections. 
This will be useful later on. 

Let us now come back to the calculation of statistical-mechanical averages. We will organize 
their calculation in two steps: first we will perform a low-temperature expansion (which is nothing 
but the Sommerfeld expansion used in the theory of free Fermi gases), expressing the finite 
temperature average in terms of a zero-temperature average. Then, we will evaluate this zero- 
temperature average by using the Wigner-Kirkwood expansion. The reason we use the low- 
temperature expansion is that, in the thermodynamic system relevant to ABJM theory, large N 
means large \i and large E, and this is equivalent to large /3, i.e. low temperature. 

We will focus on the one-body average appearing in (3.24), and we will restrict now to Fermi 
systems (i.e. we set rj = —1). We first recall that the Sommerfeld expansion expresses any 
integral of the form 

where g(E) is an arbitrary C°° function, as a power series in the temperature, 

1 = Io 9(E)dE + ^ W ( 2 " 2^) C(^ (an_1) (M)- (3-41) 
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It is easy to see that this can be written as the operator expansion, 



1 7r<9 M 

esc 



(^) e(„ - H) . (3.42) 



e /3(H- M ) + i p 
Using now (3.28), we express the average (3.24), for /? = 1, 77 = —1, as 



l(O) G0 = Tr (-^- ) = 7r^csc(7r<9> Gu), (3.43) 



where 



(3.44) 



In writing this we have used, in the first line, the fact that the star product drops out of the trace 
when only two operators are involved [32], and in the second line we used the Wigner-Kirkwood 
expansion of the distribution operator. These two expressions, (3.43) and (3.44), will be our 
basic tools to calculate vevs of Wilson loops in the Fermi gas approach to ABJM theory. 

3.2 Review of the Fermi gas approach 

We will now review the Fermi gas approach to J\f = 3 Chern-Simons-matter theories, developed 
in [46]. We will consider the generalization of ABJM theory given by necklace quivers with r 
nodes [34, 36], and with fundamental matter in each node. These theories have a gauge group 

U(N) kl x U(N) k2 x ■ ■ • U(N) kr (3.45) 

and each node will be labelled with the letter a = l, - - ,r. There are bifundamental chiral 
superfields A aa+ \, B aa -\ connecting adjacent nodes, and in addition we will suppose that there 
are Nf a matter superfields (Q a ,Q a ) in each node, in the fundamental representation. We will 
write 

k a = n a k, (3.46) 

and we will assume that 

r 

J> a = 0. (3-47) 

a=l 

According to the general localization computation in [37], the matrix model computing the 
S 3 partition function of a necklace quiver is given by 

1 f HA • ex pfe A Ll r Ui<i [2sinhf Aa -'' Aaj N )l 2 

zw ^/e^(d^H n:;L(UJj - ™ 

The building block of the integrand in (3.48) is the following iV-dimensional kernel, associated 
to an edge connecting the nodes a and b: 



i=l 



n«2coah(^) 
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Here, 



171 



U a (X) = -^X' + N fa log[2cosh 



2k 



(3.50) 



and will be interpreted as a one-body potential for a Fermi gas with N particles. We denoted by 
Aj the variables corresponding to the a node, and by pi those corresponding to the b node, after 
rescaling them as p, A — > p/k, X/k. 

We now want to interpret the kernel (3.49) as a matrix element 

K a b(Xi, • • • , Aat; Mi, • • • , Pn) = {Ai, • • • , Xn\ p a j, \p±, ■ ■ ■ ,pn} , (3.51) 

in terms of a non-symmetrized density matrix p a t> (i.e. a density matrix for distinguishable 
particles). We first notice that 

— {Al, ■ ■ ■ , Ajv| f>ab \pi, ■■ ■ , Pn} = -jyj ^2 ( _1 ) e(<J Vab (Al, • • • , \n;Pv(1), ■■■ , Va(N)) • (3.52) 



TV! 



We now use the Cauchy identity 



o-gSjv 



n i< ,[2sinh(^)] [2sinh(^) 



n M 2cosh(^) 



det 



" 2 cosh (^) 

E(-i) ew n 



(3.53) 



o-e5]v 



2 cosh 



H-Vcr(i) \ 
2 J 



In this equation, Sn is the permutation group of N elements, and e(a) is the signature of the 
permutation a. We obtain, 



K ab {X u --- ,X N ;m,--- ,p N ) = TTe v^detij ( \— —] 

^cosh^y 



A? 



k 



N\ 



(-l^JJe-^^n^ 



where we denoted 



t(x) = 



i=i 



(3.54) 



i=i 



2 cosh | 



(3.55) 



By comparing with (3.52), it follows that 

Pab(Ai,--- ,\ N ;m,--- ,p N ) = JJe _l/o(Ai) [Jt 



A? 



TV 



i=l 



i=l 



A« — /ij 
k~ 



(3.56) 



This factorization tells us that the iV-particle system is non-interacting, since the density matrix 
is completely factorized. By taking the Wigner transform of this expression, with 



h = 2irk, 



(3.57) 
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we see that p a b defines an iV-body Hamiltonian 



^ = e^ w , (3.58) 



where 

N 



i=i 

The one-particle Hamiltonian H$ is defined by 

e ; ff w =e -u a ( q ) ice -T( P ) (360) 

and 

T(p) = log (^2 cosh |) . (3.61) 

We can now use repeatedly the resolution of the identity (3.6) to write the matrix integral (3.48) 
as 

Z(N) = Tr(p), (3.62) 

where p is the density matrix 

P = P12P23 ■ ■ ■ Pr-lrPrl, (3.63) 

and this defines the one-particle Hamiltonian Hy/ by 

tt _W 12 _W 2 3 _f7 r ~ lr _bt1 

e -^w = e* " w ★ e* w ★ ■ ■ ■ ★ e* w * e* w (3.64) 

For necklace theories without fundamental matter it is easy to see that the total, one-particle 
Hamiltonian is given by [46] 

q-Hw — j; ^ \ ^ \ * . . . * \ C3 55") 

2coshf 2cosh™ 2cosh p ~ (ni 2 +W2)9 2 cosh p-("i+" 2 +^-i)<? ' 

4. Wilson loops in the Fermi gas approach 
4.1 Incorporating Wilson loops 

We will now show how the calculation of Wilson loops maps to the calculation of statitical- 
mechanical averages in the Fermi gas approach. For simplicity, we restrict ourselves to ABJM 
theory. The general N = 3 quiver can be obtained by a straightforward generalization. 

The ABJM quiver is defined by two nodes, with CS levels k and —k (as we mentioned before, 
and without loss of generality, we will take k, the level in the first node, to be positive). The 
one-body Hamiltonians associated to the edges are given by 

-H 12 i<7 2 1 -H 21 i? 2 1 

Let us consider a 1/6 BPS Wilson loop with winding number n for the first node. As shown in 
[37] and reviewed above, this corresponds to inserting 

N 



O n = ^ e nAl/fe (4.2) 
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in the matrix integral, after rescaling A — > X/k. The unnormalized vev can be written, in the 
language of many-body physics, as 

(O n ) = — J dAd/x{Ai • • • X N \O n pi 2 \m ■ ■ ■ Pn}{pi ■ ■ ■ mlhil^i • • • X N }, (4.3) 

If we integrate over p by using the resolution of the identity, we find 

(°n) = ^ J dX{X 1 ---X N \O n p 12 p2i\X 1 ---X N } = Tr(O n p 12 p 2 i) (4.4) 

which is the vev of the one-body operator (4.2) in an ideal Fermi gas of N particles with one-body 
Hamiltonian 

2 cosh V 2 cosh f l4 ' Dj 

Notice that this Hamiltonian is not Hermitian. This corresponds to the fact that the vev of a 
1/6 BPS operator is not real. After a canonical transformation, 

q^q+p, P^P, (4.6) 

we obtain a more convenient form, 

e -H w _ 1 . 1 _ e -u(g) -T{p) „x 
6 * - 2cosh|*2coshf " e+ (4 - 7j 

where T(p) is given in (3.61) and 

U(q) = log (2 cosh |) . (4.8) 

After this canonical transformation, the insertion of the BPS Wilson loop corresponds to consid- 
ering in the Fermi gas a one-body operator of the form 

O n = eXp (^), (4.9) 

which we have writen already in the one-particle sector. 

Since the vev of a 1/2 BPS Wilson loop can be obtained from (2.13) by considering the vev 
of the 1/6 BPS Wilson loop and its conjugate, we will focus on the analysis of the latter, and 
deduce the former from (2.13). 

4.2 Quantum Hamiltonian and Wigner Kirkwood corrections 

In the Fermi gas approach of [46], the full quantum Hamiltonian contains H corrections and it 
is not known in closed form. Its semiclassical expansion can be obtained by using the Baker- 
Campbell-Hausdorff (BCH) formula as applied to the *- product in (4.7). It was shown in [46] 
that, in the calculation of the grand potential of the system, only the first quantum correction 
is needed (up to exponentially small corrections in the chemical potential). However, in the 
calculation of the Wilson loop vev, we will need an infinite series of terms appearing in the 
semiclassical expansion, of the form 

U<- n \q) {T'(p)) n , rW(p) (U'(q)) n . (4.10) 
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The coefficients of these terms can be determined in closed form by exploiting some particular 
cases of the BCH formula (see [13] for examples of such calculations). We will now determine 
these coefficients. 

Let us consider the * product 

e?*ef, (4.11) 

where 

A = -U(q), B = -ap, (4.12) 

and a is a constant. The function U(q) is arbitrary. In this case, B acts as the derivative iahd q , 
and the commutator reads 

[A,B}± = iahU'(q). (4.13) 

Since this commutator commutes in turn with B, one can use a simpler version of BCH which 
says that (see for example [13]) 

Z = log, * ef ) = B + (4.14) 

where B]+ is to be understood as the operation of performing a ^-commutator with B (acting 
on the left), and its n-th power is obtained by doing the ^-commutator n times. The function 
appearing in this expression has the well-known expansion 

oo 2 

r ^ = £**" = i+f + ^ + - (««> 

n=0 

where 

c n = : (4-16) 

n! 

and B n are the Bernoulli numbers. Note that, due to a well-known property of the Bernoulli 
numbers, all the powers in this series are even except for the second term. We then conclude 
that 

Z=-ap-J2 c n {-\ha) n U^{q). (4.17) 

n>0 

Let us apply this to our case. If we take T(p) = p/2, only the terms of the form 

U^\q){T'{p)) n (4.18) 

survive in the BCH expansion. This corresponds to choosing a = 1/2 in the above formula. We 
conclude that these terms appear in the Hamiltonian in the form, 

Tip) + U(q) + ^(m n U {n) (q) {T'(p)) n . (4.19) 
n>l n - 

We can calculate in a similar way the terms obtained by exchange of p and q. 

In our case, where U(q) is given by (4.8), the derivatives of U(q) can be written in terms of 
poly logarithms. Indeed, one finds by direct computation that, 

^diTT""-^ (4 - 20) 
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where 

x = e' q . (4.21) 
Therefore, higher derivatives produce polylogarithms of lower order, 

U im \q) = (-l) m+1 Lh- m (-x), m>2. (4.22) 

Remark 4.1. Surprisingly, the real part of the Hamiltonian (4.19), when T'(p) = 1/2, can be 
written in a very suggestive form: 

P -P + H^, (4.23) 



B 2n I Irt 



where 

n=0 

On the other hand, the contribution to the free energy of the resolved conifold with g > 1 is 

oo „ 

^(-.f«) = E^^^"3_,»(e-.) 9 r 2 (4.25) 
It follows that the "quantum" Hamiltonian is related to the free energy as 

-Mf)-H)" n - xM (4 - 26) 

where 

9. = *, (4-27) 
is identified with the topological string coupling. 
Remark 4.2. By using (4.14) twice, one can compute 

-log ( e -^)/ 2 * e -°f *e-^)/ 2 ) = ap - E^T^^te)- (4.28) 

fc>0 ^ 

This determines the coefficients of (T 1 ' {p)) 2k U^ 2k \q) in the Hermitian Hamiltonian originally 
considered in [46]. The result 

-log ( e -°«*e- T W *e- a «) = 2aq + T(p) + ^(a^) 2fc ^-^|5 2fe |r( 2fc )(p) (4.29) 

fc>i 1 j - 

does not follow straightforwardly from (4.14), but it can be derived by using (4.14) together with 
symmetry arguments. This determines the coefficients of (U' (q)) 2k T^ 2k \p) in the Hermitian 
Hamiltonian in [46]. Both series of coefficients, (4.28) and (4.29), appear as well in the general 
expansion of the so-called symmetric BCH formula, and they can be verified up to high order 
with the results in [10]. The explicit, analytic expressions (4.28) and (4.29) for these coefficients 
do not seem to have appeared before in the literature. 
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We will now compute all the Wigner-Kirkwood corrections for the simplified Hamiltonian 
considered above, which is obtained from the equation 



(4.30) 



More precisely, we want to compute the generating functional of Wigner-Kirkwood corrections 
obtained by considering the Wigner transform of the canonical density matrix 

e~ tHw , (4.31) 

as explained in (3.39). To calculate (4.31), we use the following trick, inspired by similar calcu- 
lations in [13]. Let us suppose that we can write (4.31) as 

e" tH w = e - tG{q Ke- tap (4.32) 

by using the BCH formula. If this is the case, the ^-product can be easily evaluated to obtain, 

e"^w = e -tG(q) e ^J p e -tap = exp (^_ tap _ te t 9 £(g)) _ (4 33) 

where we have denoted 

£ = -iah. (4.34) 
To obtain the explicit form of G(q), we use the BCH formula (4.14) to find, 

e;' G(<?) *e-^ = exp, l_ tap _ t J2c m (tO m G^\q) \ . (4.35) 

y m>0 J 

By construction, this equals e* tHw . On the other hand, we know that 

H w = ap + ^^eU {n \q). (4.36) 



n>0 



It follows that, 



and we find 



We conclude that 



U<> /•«'/) , ^' j'un. (4.37) 



c 



l_e-€*a" w 1-e-^ 

/ e %d_ e -%d \ 

" tff w = exp -tap - - T -— iS -U(q) . (4.39) 



The second term in the exponent can be computed by using, for example, the definition of 
Bernoulli polynomials, 

T3 ^ = E^)(-ir^- (4.40) 

n>0 

One finds, 

- S °-°1fu(i) = E ""■w-^-'-'ffl/-')-^,.),,). (4 . 41) 

t l-e-Z 9 w ^ t (m + l)r w v ' 

m>0 V ' 
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Since Go = 1> Q\ = 0, the exponent in (4.39) should be of the form 

-tH w + 0(t 2 ). 



Indeed, by using 



one sees that 



B m +i(t/2) - B m+ i(-t/2) 
m + 1 



B m t + 0(t 3 



1 _ 



(4.42) 
(4.43) 
(4.44) 



The expression (4.39) generates all the functions Q r by expanding in t, and one can verify the 
results, at the very first orders, against the explicit expression in terms of ★-products given in 
(3.31). 

We also notice that the operator appearing in (4.39) can be written as 



i± d -tig 

e 2 " — e 2 " 
1 - e-59 



U(q) 



1 



1 - e-S 9 



U q + 



17 



2 



(4.45) 



One has to be however careful, since the expansion of the differential operator in the denominator 
leads to 



E 



Bt(-l) 



(4.46) 



and one has to define properly the term I = 0, which involves 1 . Comparison with the 
expansion in terms of Bernoulli numbers shows that this term stands for, 

± r [,( !+ |) -u( q - §)] -tg^i (4.4V, 
This can be also written in terms of an integral, 

2<? 



1 /■« 



17 U + 



17 u 



oo ^ 

+ g(27TT)! 



25 



(4.48) 



U {29) (A), 



where A is an appropriate reference point. 

As we will see in the next subsection, we need the expression of the canonical density matrix 
for the value 

t=f, (4.49) 

where n is the winding of the Wilson loop operator. This can be evaluated in principle with 
(4.39) and (4.46). However, the calculation is rathe delicate, since the shift by 



fi/2 = -nvri 



(4.50) 
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implies that we are resumming the series of semiclassical corrections beyond its radius of conver- 
gence, and a regularization is needed. We will proceed as follows. First, we notice that, in the 
"polygonal" limit \q\ — > oo, 

tf(9)«y + o(e~M), (4.51) 
In this limit, the second term in the exponent of (4.39) is given, for q / 0, by 

-t\q\ - |sgn(g) (4.52) 

since only U (q) and its first derivative survive. Therefore, we want to calculate the correction 
to this polygonal limit for the value of t given by (4.49). In fact, this correction is given by a 
distribution supported at q = 0, which can be obtained as follows. The first term of (4.46) can 
be computed by writing, for q > 0, 

U(q) = l + U(q), U(q) = log(l + e" 9 ). (4.53) 

We have to calculate the sum appearing in the r.h.s. of (4.47), which we write as 

1 / <<\ 2 9 



+E<^ (!)>%>■ (««> 



Since this function, as well as all its derivatives, vanish at infinity, we take A = oo as a reference 
point. For the particular value (4.50), we obtain from (4.48) 

since the integrand vanishes. We conclude that the term with i = in (4.46) is given by 

§, q>0. (4.56) 
A similar reasoning for q < shows that the first term in (4.46), for the value (4.49) of t, is 

n7ri i i /, 
-- £-M- ( 4 -57) 

The second term in (4.46) involves acting with a derivative on this first term. Equivalently, it 
can be computed as the monodromy of U(q), 

U ( Q + 7T ) - u ( Q ~ 7T J = -nvrisgn(g). (4.58) 



We then see from (4.46) that the full series of corrections involves the distribution S(q) defined 



as 



= E ^^V 1 ^ 1 sgn(g) = r -i-^ sgn (g) = M + l sgn(g) + O (0. (4.59) 
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To calculate S(q), we take a derivative w.r.t. g, and we multiply both sides by 1 — e We 
obtain the equation 

r(q)-T(q-0=S(q) + S(-q), (4.60) 

where 

T(q)=S'(q). (4.61) 

The Fourier transform of (4.60) gives 
We now set 

£ = -itf, = vrfc, (4.63) 

and solve for T(q) by doing an inverse Fourier transform. This transform is in principle ill-defined 
due to the pole at u = 0, but we can regularize it in a standard way by taking a principal value 
at the origin (or an extra derivative w.r.t. q). We obtain in this way 

T(q) = / dw e~ lujq ~. qt = - coth ( ^ ) . (4.64) 

W 2irJ sinh(^f) $ \$J V ' 



We now integrate w.r.t. q to obtain 5(g). The result is, after fixing the appropriate value for the 
integration constant, 

5(g) = i + ^log(2sinh(^)). (4.65) 

To see that this is a natural regularization, and to fix the integration constant, we note that for 
q > this can be written as 

S(q) = | + I + i log (l - e" 2 ^) , (4.66) 

while for q < we find, 

4 2 7T 7T V / 

i.e. 

5(g) = ^ + ^sgn(g) + i log (l - e 2 ^!/* 9 ) , (4.68) 
so that, for q ^ 0, and £ small, we find, 

5(<?) ~ ^ + ^sgn(g), (4.69) 

which is consistent with (4.59) and also gives the polygonal limit we need, cf. (4.52). 

It might be surprising that an infinite sum of distributions (4.59) can be resummed to a 
smooth function of q. However, this is standard in the context of the semiclassical approximation 
of Wigner functions, and is performed by means of Fourier transforms, as we have just done [55]. 
For example, the ground state of a harmonic oscillator in the Wigner formulation involves the 
Gaussian 

fw(q,p) = f(q)f(p), /(?) = ^- q2,h - (4-70) 
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Figure 2: The regions in the Fermi surface. 



But 



-Hlo 2 /4 



2vr 



. oo 



(-!)%* 



which has inverse Fourier transform 



E 

£=0 



2tt ^ 4^! 

£=o 



-co 



.2f 



(4.71) 



(4.72) 



In the classical limit h — > we have a localized particle at the origin, and the h corrections give 
an infinite sum of distributions which however can be obtained from the smooth Gaussian in 
(4.70). 

We conclude that, for the simplified Hamiltonian (4.30), the canonical density matrix with 
t given in (4.49) is 



2n tj 



exp 



n mri ( ( 

p H — n log I 2 sinh I 



k c ' 2 

at least with the natural regularization procedure explained above 



(4.73) 



4.3 Integrating over the Fermi surface 

We are now ready to calculate the vev of the 1/6 BPS Wilson loop with winding number n in 
the Fermi gas approach. The corresponding one-body operator is given in (4.9). The first step 
is then to calculate (3.44) for this operator, i.e. 

= 7 - ^ w)e ^ + E H d^-i y ii 6 ^- H ^re * • (4.74) 

r>l 
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Notice that, since the Hamiltonian is complex, the Fermi surface 



H w (q,p) = » (4.75) 

is in principle a surface in complex space. However, up to exponentially small corrections, we 
can recover a real Hamiltonian by a Wick rotation of the Planck constant, h — > —ih, so that ih 
is real. After doing this, the integration process is perfectly well defined, and we can rotate back 
at the end of the calculation. Equivalently, it can be easily seen from our computations that this 
involves integrating over appropriate paths in the complexified phase space. 

The first integral is over the region enclosed by the Fermi surface. However, by integrating 
w.r.t. p or q one can reduce the integral to a boundary integral over the Fermi surface, plus a 
"bulk" contribution which is easy to calculate. As in [46], we will divide the boundary of the 
Fermi surface in appropriate regions. The quantum Hamiltonian reads, 

H W = T{p) + U(q) + ^U'(q)T'(p) + ■■■ (4.76) 
where the corrections are exponentially small. The point in the Fermi surface with p coordinate 

P* = n + - (4.77) 



has a q coordinate of the form 



q* = li+^ + 0{e-»). (4.78) 



It is easy to see that the leading contribution to the Wilson loop is obtained by subtracting the 
contribution of the bulk region 

-p*<p<p*, -q*<q<q* (4.79) 

to the contribution of the boundary (of course, in writing this inequalities, we assume that we 
have performed a Wick rotation and that ih is real). But, if we restrict ourselves to terms which 
are proportional to exp(2n/x/&), the only contribution comes from the boundary shown in red in 
Fig. 2. This region can be divided in turn in two regions: a region where 

p > p*, <q<q*, (4.80) 

and the region obtained by exchanging p and q, 

q>q*, ~P*<P<P*- (4-81) 

They give the same contribution, so we will restrict ourselves to the first region and then multiply 
the result by two. Along the curve bounding the region (4.80) we can neglect exponentially small 
terms in p, i.e. we can assume that T(p) = p/2. We can then write 

p(ji,q) = 2 f i + (2H w -p), (4.82) 

where 

2H w -p = U(q) + jU'(q) + --- (4.83) 
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only depends on q and it has been computed in (4.19), with T(p) = p/2. We want to calculate 
the first term in (4.74), 

f dqdp nig+p) k f dq ? / n P ^, q ) \ 

J—e u ei/i -H w ) = -J—eu{e * - l) (4.84) 

and we restrict to terms which are proportional to exp(2n^//c), so we keep only the first term. 
After plugging in the value of p(fi, q), we find 



k 2 
e 



2itnh 



— / dge^^e"^ w . (4.85) 

J — q* 



Notice that the p dependence in this and similar expressions cancels trivially. It is easy to see that 
all h corrections to the Hamiltonian contribute to this integral, even if we neglect exponentially 
small corrections. 

Let us now consider the Wigner-Kirkwood corrections to (4.74) along the curve bounding 
the region (4.80). By writing 

S(fi - H w (q,p)) = dH l {qp) M P ~ ?)) = 25(p - P(M, ?)) (4-86) 



we obtain 



(4.87) 



This combines with (4.85) to produce, 

(4.88) 

Using (4.73) we find that this integral equals 

k „ 2n M fi* e n i/ k 

1 e ^ / dg , . - , q ,, w , (4.89) 



^nh J_ q , y (2sinh(f)) 

where 

g»»/i + ^-- ( 4 - 90 ) 

There is a singularity of the integrand for g = 0. However, as we will see, this can be avoided in 
a natural way. Also notice that, as g* — > oo, the integral diverges due to the upper integration 
limit, but it is not divergent when we send the lower integration limit to infinity. In fact, doing 
this only introduces exponentially small corrections (which we are neglecting anyway), and up 
to these corrections we can just compute, 

/q* & nq/k 
. dg . w . (4.91) 

-oo-^ ( 2sinh (f)) 
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To calculate this integral, we make the following change of variables 

u = e q/k , (4.92) 

so that 

In= / du- (4.93) 

It is actually simpler to calculate the generating functional, 

where 

and we have again neglected exponentially small corrections. We now take into account that 

log(l - z) 



= Y J H n z n , (4.96) 



1 ~ z 

n=l 

where H n are harmonic numbers, to obtain 

M n = »-'^- k -H n _ x . (4.97) 

Notice that the integrand above has poles at u 2 = (1 — z)^ 1 , and we have chosen an integration 
contour in the complex u-plane which avoids these poles. This is natural since the upper limit 
of integration, u*, is in fact complex. 
Putting all together, we obtain 

e k i I fx - --ffn-i • (4.98) 



2irnh 



As we explained above, there is an identical contribution from the region obtained by exchanging 
p o q. Finally, one has to subtract the contribution from the bulk region, which gives 

f q * f P * dqdp "(p+g) k 2 2n M ,inH i n k 2 2n M 

" U ^ e ^~ = + • • • = "« e " + • • • (499) 

where the dots denote subleading exponentially small corrections. Therefore, up to these correc- 
tions, we find 

no n (p) « ^e^i" (V - ^ - kH r ) . (4.100) 



2imh V 2 

As we will see in a moment, this is in precise agreement with the result obtained in the 't Hooft 
expansion at genus zero. 

According to (3.43), in order to find the full statistical- mechanical average, we just have to 
take into account the finite temperature corrections encoded in the Sommerfeld expansion. We 
then find, 

l{O n ) GC = ndncaeiirdJnoM, ( 4 - 101 ) 
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with the value obtained in (4.100), which we will write as 

n oM ~ (A(k)fi + B(k))e^ 
Here A(k) and B(k) are given by 

k 



Putting things together, we find 

^(^) GC = 2 | n cs C 2 | n 



B(k) = - 



k 



:n+l 



4ir 2 n 



TT 



(fi+^-ncot 27 ^)A(k) + B{k) 

ATI K 



e k 



where S is the grand potential calculated in [46] 

„ / 2/i 3 [i \ik\ 

" = eXP fe + 3fe + 24j' 



(4.102) 



(4.103) 



(4.104) 



(4.105) 



up to exponentially small corrections and an overall, /U-independent constant. (4.104) gives then 
the exact grand canonical correlator at all k, up to exponentially small corrections in //. To get 
the original normalized Wilson loop correlator, we have to perform the inverse transformation 



where Z denotes the partition function of the theory, which is itself given by 

Z(N) = ^-J' d/xe-^-Ou). 
We recall that the Airy function has the integral representation 



Ai(z) 



1 
2^i 



dt exp 



t 6 



zt 



(4.106) 



(4.107) 



(4.108) 



where C is a contour in the complex plane from e 17r / 3 cx) to e 17r / 3 oo. Therefore, 

k 



Z(N) oc Ai C- 1/? '[n 



24 3k 



-) 

Ik) 



(4.109) 



which is a result first derived in [23] for ABJM theory and then rederived in the Fermi gas 
approach in [46] for a class of N = 3 theories. Now, due to the exponential form of (4.104), the 
integral (4.106) can be written in terms of quotients of Airy functions, 
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where the prime denotes the derivative of the Airy function, and 
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The functions A\(k) and A 2 (k) are defined as 
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Once the answer for the 1/6 BPS Wilson loop correlator is found, we can obtain the expectation 
value of the 1/2 BPS Wilson loop via (2.13) 
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Notice that the Airy functions in the denominators of (4.110) and (4.114) come from the partition 
function [23, 46]. 

We should emphasize that (4.110) and (4.114) are exact results at all orders in the 1/N 
expansion, up to exponentially small corrections. We will now extract from it some results on 
the 't Hooft genus expansion and test it with known results at low genus. 

4.4 Genus expansion 

In order to extract the 't Hooft expansion of the Wilson loop correlator, we have to expand 
(4.110) in powers of 1/k. Since we are working in the 1/N expansion and k is generic, the results 
we will obtain are valid in the strong 't Hooft coupling regime 



A > 1 



(4.115) 



The 't Hooft expansion of the Wilson loop vevs from the ABJM matrix model has been reviewed 
and extended in section 2. Therefore, we can compare the genus expansions we obtain from 
(4.110) with known results. We will do the expansions explicitly for genus zero, genus one, and 
genus two. In appendix B, we will summarize the results for few more higher genus expansions. 
In the literature, the normalized expectation values of Wilson loops are further normalized by 
the dimension of the representation of the Wilson loop. Since we are considering Wilson loops 
in the fundamental representation, we will divide (4.110) by an extra factor of N. 



— Genus zero 

To test the agreement between the results of ABJM matrix model and the Fermi gas ap- 
proach, it will be more convenient to expand the Airy functions in (4.110) in terms of the 
/•c-variable (2.53) where only positive powers of k are relevant at strong coupling. For genus zero 
we find, 

< w 2 / *>'— ^(""•"-T-"-)- (4116) 

This agrees with the result (2.51) obtained with standard matrix model techniques. 
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The strong coupling expansion of this result can be obtained by expanding the Airy functions 
in (4.110) in terms of the 't Hooft coupling A in the regime (4.115), 
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Once we have obtained the result for the expectation value of the 1/6 BPS Wilson loop, the 
result for the 1/2 BPS Wilson loop follows from (2.13) by incorporating the result of the other 
node of the ABJM quiver gauge theory. 

— Genus one 

As the next step of our checks, we would like to compare the results of the Fermi gas approach 
and the ABJM matrix model at genus one. Expanding (4.110) in terms of k and picking the 
coefficient of k~ 2 , we find the following expression 
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Using (2.13) we obtain for the vev of the 1/2 BPS Wilson loop at genus one, 



(4.118) 
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The genus one result for the 1/6 Wilson loop correlator with winding one was first found in [15] 
by analyzing the ABJM matrix model. If we set n = 1 in (4.118), we find precise agreement 
between our results and those of [15]. We can also easily compute the genus one result for the 
1 /2 BPS Wilson loop correlator with arbitrary winding from the ABJM matrix model, following 
[15], and the result is in agreement with the general expression (4.119) 

Expanding the Airy functions in (4.110) directly in terms of the 't Hooft coupling A in the 
region (4.115), we find the following expansion for the 1/6 BPS Wilson loop expectation value 
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Similar to the genus zero result, the expansion of the 1/2 BPS Wilson loop correlator with 
winding n is automatically obtained by applying (2.13). 



— Genus two 

As our last check, we consider ABJM Wilson loop correlators at genus two. Expanding 
(4.110) in terms of re, we have 
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For n = 1, the above expression specializes to 
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Needless to say, the matrix model computations at genus two are very involved. This expression 
agrees with the result of the matrix model calculation at g = 2 and strong coupling. The 't Hooft 
expansion at strong coupling at genus two is found from (4.110), 
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where 1/2 BPS Wilson loop result at genus two is immediately obtained by applying (2.13) to 
(4.123), as in the previous cases. 

Since we have the exact result (up to exponentially small corrections) for the Wilson loop 
correlator (4.110), we can extract the leading and next to the leading terms of the 't Hooft 
expansion at arbitrary genus and strong coupling, as it was done in [14] for the 1/2 BPS Wilson 
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loop of = 4 super Yang-Mills theory. For the 1/6 BPS Wilson loop correlator with arbitrary- 
winding n, we find 
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Notice that the above expressions are normalized to k~ 29 . If one desires to normalize to gl 9 ', one 
needs to divide by an extra factor of (27r) 29 (— I) 9 . 

Using (2.13), the leading and next to the leading terms of the 1/2 BPS Wilson loop correlator 
are found 
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It turns out that the 1/2 Wilson loop correlator does not involve c g coefficients. At every genus, 
the ratio of the leading terms of the 1/6 and 1/2 Wilson loop expectation values is given by 
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This ratio was first found in [45] at genus zero for the trivial winding, and (4.128) generalizes 
this result for arbitrary genus and winding. 

There are two interesting properties of the expression (4.110) and (4.114) which are worth 
noticing. First of all, both expressions are singular for 



k = 1, 2. 



(4.129) 



Indeed, for these values of k, the semiclassical expression for the vev (3.24), which is simply 
an integral over phase space, is not convergent, and our final result is reflecting this through a 
pole in the esc function. We conclude that our expression is not valid for these special values 
of k. This is not that surprising, since precisely for those values, ABJM theory has enhanced 
supersymmetry to N = 8 [27]. On the other hand, and as we have just tested, the Fermi gas 
result resums the genus expansion for k large, therefore the value k = 2 sets the convergence 
radius for this expansion and the resulting function has the above singularities. It would be very 
interesting to understand what happens when k = 1, 2. 
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Second, we recall that the vevs of 1/2 BPS Wilson loops in ABJM theory can be related 
to open topological string amplitudes in local P 1 x P 1 [45, 15]. The expression (4.114) can be 
interpreted as saying that the unnormalized Wilson loop vev is given by the Fourier transform 
of a disk amplitude at large radius, 



where n is the multicovering degree. As pointed out in [46], the Fourier transform expression 
(4.107) for the canonical partition function can be interpreted as the change of symplectic frame 
from the large radius topological string partition function, to the orbifold partition function [1]. 
The result (4.114) indicates that a similar result should be valid for the open sector, namely, 
that changes of frame in the open sector are also implemented by Fourier transforms of the open 
string partition function. 

We can now identify the overall inverse sine appearing in this formula: it is nothing but the 
well-known all-genus bubbling factor for a disk in topological string theory. This was first found 
in [50] by using large N duality with Chern-Simons theory and derived in [38] with localization 
techniques. Amusingly, we have re-derived this factor here by using Sommerfeld's expansion for 
Fermi gases at low temperature. This bubbling factor resums the genus expansion, as in the 
Gopakumar-Vafa representation of the closed topological string free energy [24]. We see however 
that this resummation leads to singularities at k = 1,2, and this was already observed in the 
closed string sector in an attempt to resum worldsheet instantons a la Gopakumar-Vafa [47]. 

5. Conclusions and prospects for future work 

In this paper we have used the Fermi gas approach developed in [46] to compute vevs of Wilson 
loop observables in the ABJM matrix model of [37]. This approach is based on reformulating 
the matrix integral as the partition function of an ideal Fermi gas in an exterior potential, and 
on a semiclassical evaluation of the resulting quantities. The calculation of Wilson loop vevs is, 
however, more difficult than the one of the canonical partition function made in [46] , since one has 
to resum an infinite number of corrections in h. We have seen however that such a resummation 
is feasible and we have obtained expressions valid at all orders in the genus expansion and for 
strong 't Hooft coupling. Equivalently, the expressions we have obtained are full M-theoretic, 
since they are valid for finite k and large N. 

Clearly, this work can be generalized in many different ways. Already in ABJM theory, 
it would be interesting to extend our calculation to higher representations. The Wilson loop 
operator for a representation with £ boxes is an ^-body operator in the Fermi gas. Let us 
consider for example a 1/6 BPS Wilson loop in the antisymmetric representation. We can write 
the corresponding matrix model operator as 




(4.130) 




(5.1) 
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so it is clearly a two-body operator. We then have, schematically, 




GC 



1 




2 



(5.2) 




- 41 - 



which is a sum of "direct" (Hartree) and "exchange" (Fock) terms. The first term factorizes into 
the product of two one-body operator vevs, like the ones computed in this paper, and we are 
left with the calculation of the exchange term. This can in principle be computed as well with 
semi-classical techniques (see [48] for a closely related example) , therefore it would be interesting 
to develop these techniques further. 

Another obvious generalization of this work would be to consider Wilson loop vevs in the 
M = 3 Chern-Simons-matter theories analyzed from the Fermi gas point of view in [46]. In 
these theories it is general difficult to obtain results by using the traditional tools of matrix 
models in the 't Hooft expansion, therefore the procedure developed in this paper is probably 
the simplest one to go beyond the large N limit. However, a detailed analysis will require 
computing quantum corrections to the Hamiltonian, Wigner-Kirkwood corrections, and finally 
find an appropriate regularization of the resummed semiclassical expansion. In fact, it would 
be interesting to understand in more detail the regularization procedure developed in section 4. 
In ABJM, this procedure could guided from the results in the 't Hooft expansion, but for more 
complicated M = 3 theories we might need a better understanding of this issue. 

Our final result for 1/2 BPS Wilson loops in terms of Airy functions(4.114), gives an analytic 
result for finite k ^ 1,2 and large N, at all orders in l/N, and up to exponentially small 
corrections. Such an expression is very well suited for the type of numerical testing performed 
in [29]. In that paper, numerical calculations provided a beautiful verification of the analytic 
formulae proposed in [15, 23, 46], and helped in clarifying certain aspects of the analytic results 
(like for example the nature of the function A{k) introduced in [46]). Such a numerical test would 
be also very useful in understanding what happens when k = 1, 2, where our formula displays a 
singular behavior. 

Finally, we pointed out that (4.114) has a natural interpretation as a generalization of the 
Fourier transform of [1] to the open sector. It would be clearly very interesting to understand this 
better, and more generally, to develop techniques to compute topological open string amplitudes 
at higher genus in an efficient way. 
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A. 1/6 BPS Wilson loops at arbitrary winding number 

In this appendix, we present the details of the matrix model computation which led to (2.49). 
Our starting point is the integral (2.48). This integral can be explicitly evaluated in terms of 
elliptic functions [8] 




(A.l) 
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where the functions Vj are defined recursively, in terms of elliptic integrals, as follows 

V = K(k), 
V 1 = U(a 2 ,k), 

V 2 = ^ \o?E{k) + (k 2 - a 2 )K(k) + (2a 2 k 2 + 2a 2 -a 4 - 3k 2 )U(a 2 , k)] , 
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V " = 2 (m + 2)(l-^)(fc2-^) [(2m + l)k2Vm + (2m + 1)(a2fc2 + a2 ~ 3fc2)Kl+1 
+(2m + 3) (a 4 - 2a 2 A: 2 - 2a 2 + 3k 2 )V m+2 ] . 

(A.2) 

In the above expressions, the moduli of the elliptic functions are defined by 

jfc , = (.'-l)tf»-l) q2 = 1z4. (A .3) 
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To understand the strong coupling behavior of the integral (2.48), we need to expand the above 
functions in the large k regime. First, we notice that 
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We can easily find the solution to the above recursion relation. We have 
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Using the above solution (A. 6), we then compute the integral (2.48) in regime of large re 

^-iV-^E^^ + f + |log/c)(-l)^Q)+2Iog«j . (A.7) 

In order to perform the sum on the harmonic numbers in the above expression, we use the 
following integral representation of harmonic numbers 

H S - X = / — — dx . (A.8) 
Jo 1 ~~ x 

Using the above representation, we then obtain 

= |^(l-(l- x)- 1 ) = H n _, . (A.9) 
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The sum on the rest of the terms in (A. 7) is easy to perform 

E(-1) J (j) (f + ^og«) +21og^ = -j + \log* . (A.10) 
Putting things together, we therefore conclude 
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which is the sought-for result (2.49). 

B. Results at g — 3 and <? = 4 

It is evident from (4.110) that extracting higher genus contributions to the expectation values of 
the ABJM Wilson loops is more economical than other existing approaches, such as the matrix 
model approach. To demonstrate this, we summarize the result of the 't Hooft expansion for 
g = 3 and g = 4 in this appendix. It will not be necessary to find the expansions in terms of k, 
as there is no matrix model computation that we can compare against it. It will be sufficient to 
directly expand the Airy functions in (4.110) in terms of the 't Hooft coupling at strong coupling 
regime. 

For g = 3, we obtain 
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while the g = 4 contribution is given by 
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Applying (2.13), we can immediately find the result for the 1/2 BPS Wilson loop correlator at 
5 = 3 and g = 4. 
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